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BACKGROUND 


As  part  of  the  Submarine  Technology  Program,  the  Defense  Advanced  Research  Projects 
Agency  (DARPA)  recently  sponsored  a  Low-Frequency  Structural  Acoustics  Benchmark 
Exercise.  The  purpose  of  this  exercise  was  to  test  and  validate  several  major  computational 
codes  that  have  been  developed  to  solve  acoustic  scattering  problems  of  elastic  objects  in  a 
fluid.  All  of  the  computations  of  scattering  were  done  on  a  Qay  YMP  at  the  Los  Alamos 
National  Laboratory  (LANL).  The  result  of  each  problem  was  a  large  digital  data  set,  which 
was  analyzed  and  displayed  off-line  on  a  workstation  using  specially  developed  visualization 
techniques.  The  Benchmark  Exercise  began  in  May  1991,  and  a  comprehensive  report  on 
the  results  will  be  issued  separately  by  LANL. 

The  purpose  of  this  report  is  to  describe  some  of  the  visualization  techniques  and 
procedures  that  were  developed  to  review,  compare,  and  analyze  the  large  amount  of 
computational  data  generated  in  the  Benchmark  Exercise.  It  was  felt  that  the  visualization 
approach  and  techniques  should  be  described  in  an  independent  document  because  the 
techniques  are  generally  useful  in  representing  scattering  or  target  strength  functions 
generated  by  any  means,  and  the  approach  is  applicable  to  any  computational  or  experimen¬ 
tal  problem  in  which  there  is  a  need  to  understand  large  multidimensional  data  sets. 
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THE  NEED  FOR  VISUAUZATION 


In  the  Benchmark  Exercise,  a  series  of  nine  problems  was  developed  to  test  the  accuracy 
and  computational  efficiency  of  six  different  structural  acoustics  codes.  In  each  of  the 
problems,  the  forcing  function  was  a  steady-state  plane  wave  at  a  set  of  densely  spaced 
frequencies  incident  upon  the  scatterer  from  one  of  several  different  directions.  In  all  the 
problems,  it  was  required  that  the  code  generate  the  complex  numbers  representing  the 
complete  three-dimensional  far-freld  scattering  frinction.  In  some  of  the  problems,  it  was 
also  required  that  the  code  generate  the  complex  surface  pressure  and  normal  velocity  on  a 
dense  grid  over  the  complete  surface  of  the  scatterer.  The  problems  and  results  were  specified 
in  nondimensional  terms. 

A  far-field  data  set  typically  consisted  of  all  the  information  needed  to  reconstruct  the 
normalized  target  strength  or  complex  scattering  function  at  a  set  of  frequencies  over  all 
space.  Each  set  was  given  a  name  associated  with  the  problem  and  the  code  that  generated 
it,  and  an  extension  of .  f  f  #,  where  the  #  was  a  single  digit  identifying  the  incidence  angle. 
Each  set  was  a  file  containing  a  three-dimensional  array  of  complex  coefficients  in  a 
circumferential  harmonic  decomposition  of  the  three-dimensional  field,  as  follows: 


N-l  N-l 

P(fp)  •  ^  anCos(nq>)  +  ^  b»sw(n«p), 
nmO  rimO 


where  N  is  the  number  of  harmonics  and  bo  -  0  .The  dimensions  of  the  array  are  the  number 
of  frequencies,  the  number  of  azimuthal  harmonics,  and  the  number  of  bistatic  observation 
angles.  In  a  typical  problem,  the  computation  was  made  at  331  frequencies,  73  observation 
angles,  and  3  incidence  angles,  using  7  harmonics.  Each  complex  number  required  a  16-bit 
(IEEE)  binary  representation,  so  the  far-field  output  from  one  problem  typically  comprised 
three  1.4-megabyte  files. 

Similarly,  a  surface  pressure  (* .  sp#)  or  normal  velocity  (* .  nv#)  data  set  was  also  a 
three-dimensional  array  of  complex  numbers,  whose  dimensions  were  the  number  of 
frequencies  by  the  number  of  harmonics  by  the  number  of  grid  points  describing  the 
generator  line  of  an  axially  symmetric  object.  The  number  of  grid  points  varied  with  the 
code,  but  was  typically  of  order  250  in  a  large  problem,  resulting  in  surface  files  of 
approximately  5  megabytes  each. 

Some  of  the  problems  were  run  more  than  once  for  each  structural  acoustic  code,  so  the 
total  amount  of  computational  data  generated  and  reviewed  in  the  Benchmark  Exercise 
approached  5  gigabytes.  Clearly,  an  efficient  and  powerful  graphical  postprocessing  method 
was  required. 
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APPROACH 


In  addition  to  the  sheer  magnitude  of  the  problem,  there  were  many  other  challenges 
inherent  in  trying  to  display,  compare,  and  understand  the  results.  New  formats  had  to  be 
conceived  and  implemented  to  display  or  compare  large  amounts  of  data  at  one  time.  Hie 
judicious  use  of  color  and  animation  to  make  the  results  more  understandable  was  an¬ 
ticipated.  In  addition  to  rapid  and  convenient  software  development,  the  need  for  efficient 
repetition  of  tasks  for  many  different  data  sets  was  recognized.  Reproduction  of  workstation 
displays  on  paper  or  videotape  was  also  important. 

With  these  criteria  in  mind,  the  command  language  version  of  PV-WAVE*  software 
was  chosen  as  the  means  for  accomplishing  the  visualization  tasks  on  available  workstations 
(Sun  and  Silicon  Graphics).  The  key  features  of  this  software  that  were  of  importance  to  the 
task  at  hand  were  as  follows: 

•  Availability  of  a  large  number  of  basic  mathematical  and  graphical 
operations  as  simple  callable  procedures  and  functions 

•  Availability  of  a  high-level  programming  language,  with  ail  the  usual 
logical  constructs 

•  Ability  to  run  the  software  interactively,  retaining  data  and  intermediate 
variables  in  memory 

•  Ability  to  generate  new  user-defined  procedures  and  functions  for  efficient 
development  of  a  family  of  visualization  algorithms 

•  Ability  to  develop  compiled  user-friendly  routines  for  production  use  by 
junior  programmers 

•  Speed  and  efficiency  of  operation  on  a  workstation,  especially  with  matrix 
or  array  operations 

The  example  programs  in  this  document  are  all  written  in  the  PV-WAVE  command 
language,  but  should  be  understandable  by  anyone  familiar  with  a  high-level  structured 
programming  language. 


*  PV-WAVE  is  the  trademark  for  a  software  product  of  Precision  Visuals,  Inc,  Boulder,  CX). 
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ARCHITECTURE  AND  DATA  FLOW 


An  overall  guide  to  the  approach  used  is  described  in  figure  1.  At  the  top  of  the  flowchart 
are  shown  the  IEEE  binary  data  files  written  to  disk  by  the  structural  code  that  was  executing 
on  the  Cray  YMP.  All  binary  files  began  with  a  header  in  which  problem  identification  and 
parameter  descriptions  as  used  at  the  time  of  computation  were  written.  In  addition  to  the 
far-field,  surface  pressure,  and  normal  velocity  files  described  above,  an  additional  binary 
coordinate  (* .  boo)  file  was  written  for  each  problem,  containing  the  set  of  cylindrical 
coordinates  describing  the  shape  of  the  scatterer  in  discretized  form.  These  coordinates  were 
also  the  locations  at  which  the  surface  pressure  and  normal  velocities  were  calculated. 

Each  shaded  rounded  rectangle  in  figure  1  represents  one  or  more  PV-WAVE  command 
language  routines  written  to  perform  a  specific  function.  These  routines  are  small  ASQI 
files  having  a  .  pro  extension,  created  with  any  convenient  word  processor.  On  the  top 
level  of  the  flowchart,  there  are  three  conversion  routines  (convf fb.pro,  con- 
vbco .  pro,  and  convsurb .  pro).  These  and  other  .  pro  routines  described  in  this 
section  are  listed  in  their  entirety  in  appendix  A.  The  conversion  routines  are  intended  to  be 
used  only  once  on  a  given  data  set  The  purpose  of  each  conversion  routine  is  to  read  an 
IEEE  binary  file  generated  by  a  FORTRAN  program  on  the  Cray  YMP,  and  rewrite  it  as  a 
structured  binary  file  in  a  format  optimized  for  input  to  PV-WAVE.  One  of  the  built-in 
functions  in  the  PV-WAVE  language  is  the  ASSOC  function,  whose  purpose  is  to  map  an 
array  structure  onto  a  named  file,  thereby  permitting  the  entire  file  to  be  read  efficiently  into 
memory  and  to  be  associated  with  a  structured  variable  in  one  statement.  Once  the  structured 
binary  files  (* .  f  w#  for  far-field  PV-WAVE  files,  * .  pw#  for  surface  pressure,  and  * .  vw# 
for  normal  velocity)  are  created,  the  original  IEEE  binary  files  can  be  discarded  or  at  least 
archived  because  they  are  not  needed  for  subsequent  processing  of  the  data. 

The  convbco .  pro  routine  extracts  the  coordinate  information  from  the  * .  bco  file, 
writes  it  to  a  *  .wco  file,  and  also  extracts  the  header  information  and  writes  it  to  an  ASQI 
* .  hdr  file.  This  header  file  is  the  key  to  generalizing  the  remainder  of  the  visualization 
routines  so  that  a  single  simple  * ,  pro  routine  will  work  for  any  problem,  with  any  choice 
of  modeling  variables,  independent  of  the  particular  structural  acoustics  code  that  generated 
it.  The  information  in  the  header  is  made  available  by  calling  the  procedure  rdhdr .  pro 
and  passing  the  information  as  named  variables  through  a  common  block  in  rdhdr .  pro 
and  the  routine  that  calls  it.  Calling  rdhdr .  pro  also  causes  the  header  information  to  be 
displayed  on  the  screen  ,  as  shown  in  figure  2.  The  five  types  of  PV-WAVE  files  on  the 
third  level  of  figure  1  contain  all  the  information  present  in  a  solution  to  one  of  the  benchmark 
problems. 

For  the  purpose  of  quickly  checking  the  data  set  or  for  developing  some  new  means  for 
visualizing  it,  there  are  two  workhorse  routines  called  grabf  f .  pro  and  grabsur .  pro. 
The  purpose  of  these  routines  is  to  quickly  read  (i.e.,  “grab”)  a  far-field  or  surface  variable 
data  set  and  its  associated  header,  and  make  the  information  available  in  high-speed  memory 
for  use  on  the  workstation.  Each  routine  also  displays  a  simple  representative  portion  of  the 
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data  on  the  screen.  In  the  case  of  grabff  .pro,  this  representation  is  a  set  of  three 
two-dimensional  color  images,  as  shown  in  figure  3.  One  of  these  (figure  3a)  shows  target 
strength  in  dB  encoded  in  color  as  a  function  of  polar  angle  and  nondimensional  fi-equency 
for  a  fixed  azimuthal  angle.  The  other  two  surfaces  (figures  3b  and  3c)  show  the  linear 
magnitude  of  the  scattering  function  and  its  phase  in  the  same  layout.  A  color  table  is 
included  with  each  figure  to  show  how  the  tabulated  values  are  translated  into  a  color  image. 
Horizontal  cuts  of  these  surfaces  represent  directivity  patterns  at  some  fixed  frequency,  and 
vertical  cuts  represent  spectral  variations  in  target  strength  at  some  (possibly  bistatic)  angle. 
In  the  case  of  grabsur.pro,  two  color  images  are  given  in  figure  4:  a  “Persian  rug” 
(figure  4a),  in  which  surface  pressure  or  normal  velocity  is  plotted  as  a  function  of  axial  and 
azimuthal  coordinates  for  a  fixed  frequency,  and  a  “standing  wave”  plot  (figure  4b),  in 
which  the  axial  variation  of  the  surface  quantity  is  displayed  as  a  function  of  axial  coordinate 
and  frequency  for  a  fixed  azimuthal  angle.  If  there  are  any  data  format  problems,  or  if  a 
fundamental  error  has  been  made  in  computing  the  results,  this  initial  visual  check  is  often 
sufficient  to  identify  the  problem.  Using  the  built-in  procedures  and  functions  of  PV-WAVE, 
the  data  in  memory  can  be  easily  manipulated  and  presented  for  visual  analysis  on  the 
workstation  screen  in  a  variety  of  one-,  two-,  and  three-  dimensional  formats.  This  is  a  fast 
and  convenient  way  to  explore  the  data.  In  addition,  PV-WAVE  provides  both  a  “command 
history”  and  a  “journaling”  feature  that  make  this  interactive  mode  very  efficient  for  the 
development  of  new  algorithms.  Instructions  can  be  issued,  recalled  and  modified,  and  then 
saved  to  form  the  core  of  a  new  *,pro  routine,  because  the  syntax  of  the  language  is 
virtually  identical  whether  used  interactively  or  read  and  compiled  fi-om  a  program  file. 

The  rectangle  labeled  other .  pro  is  a  generic  descriptor  on  the  chart  for  all  the  other 
compiled  routines  that  are  described  in  this  document  or  that  were  developed  during  the 
Benchmark  Exercise.  A  number  of  such  application-oriented  routines  are  listed  in  the  next 
section,  along  with  a  brief  statement  of  the  scope  or  purpose  of  the  routine.  Examples  of  the 
output  are  also  given  as  appropriate,  although  these  outputs  are  not  a  direct  representation 
of  what  would  be  displayed  on  the  screen  by  the  PV-WAVE  routine.  The  complete 
PV-WAVE  code  for  these  routines  given  in  appendix  A;  however,  it  is  not  intended  that  the 
code  could  be  used  as  is,  since  in  many  cases  it  is  dependent  on  the  particular  data  sets 
available  for  analysis. 

Finally,  in  reference  to  figure  1,  in  order  to  redirect  the  output  of  the  various  routines  to 
a  printer  or  other  hardcopy  device,  two  simple  routines  called  hardcopy .  pro  and 
closeit.pro  were  written.  These  routines  can  be  used  interactively  or  inserted  into 
compiled  visualization  code  surrounding  one  or  more  graphical  display  statements.  They 
reformat  and  redirect  the  output  to  an  encapsulated  PostScript  file  for  inclusion  in  a  document 
or  viewgraph  using  some  convenient  page-formatting  software.  For  production  purposes, 
some  routines  were  rewritten  to  make  the  hardcopy  the  normal  output,  or  to  do  both  so  that 
the  process  could  be  monitored  as  it  was  generating  the  files  for  producing  hardcopy  output. 
Videotapes  were  also  generated  using  some  routines  (especially  for  those  cases  where 
animation  in  time  was  useful  for  analysis  of  the  effects  of  changes  in  frequency  or  angle). 
The  simplest  way  to  generate  the  videotape  is  to  connect  appropriate  scan  converting  and 
recording  equipment  directly  to  the  workstation. 
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EXAMPLE  APPLICATIONS 


Selected  examples  of  visualization  and  analysis  routines  are  presented  in  this  section. 
The  first  two  examples  are  routines  used  simply  to  display  or  visualize  far-field  and  surface 
results.  The  next  two  examples  were  developed  to  aid  in  the  comparison  of  multiple  data 
sets.  The  last  two  routines  actually  process  the  results  in  such  a  way  as  to  reveal  the 
underlying  physical  mechanisms. 

f  f  3d.pro  -  The  purpose  of  this  routine  is  to  fully  use  the  graphic  display  device  to 
show  the  three-dimensional  character  of  the  far-field  scattering  function,  using  both  color 
and  distance  from  the  origin  to  indicate  target  strength  at  a  given  frequency.  This  gives  a 
more  complete  geometrical  picture  of  the  scattering  function  than  the  color  images  produced 
by  grabf  f .  pro  which  show  the  variation  in  polar  angle  only,  for  a  fixed  azimuthal  angle. 
This  visualization  techinique  is  particularly  useful  for  understanding  whether  the  results 
have  the  correct  symmetry  in  a  particular  problem.  Because  each  three-dimensional  picture 
is  for  a  fixed  frequency,  the  variation  in  frequency  is  shown  by  animating  the  display;  that 
is,  by  showing  a  sequence  of  frames  at  some  chosen  frequency  increment.  Figure  5  shows 
several  frames  in  such  an  animation. 

plot  sur « pro  -  The  purpose  of  this  routine  is  to  overlay,  by  means  of  animation,  line 
plots  of  either  surface  pressure  or  normal  velocity  as  a  function  of  axial  node  number.  This 
is  useful  to  determine  the  spatial  variation  in  relation  to  the  surface  or  internal  features  of 
the  scattering  object,  and  to  observe  and  understand  the  difference  between  pressure  and 
velocity  on  the  surface.  Figure  6  indicates  the  final  frame  of  the  animation,  which  is  built 
up  over  time  by  overlaying  many  different  velocity  distributions  in  a  selected  frequency 
range. 

envel.pro  -  The  purposes  of  this  routine  are  to  generate  a  set  of  overlays  of 
supposedly  equivalent  solutions  generated  by  different  structural  acoustic  codes,  to  allow 
the  user  to  select  a  subset  of  the  solutions  (possibly  having  rejected  one  or  more  outliers), 
and  then  to  produce  the  mean  and  an  envelope  that  bounds  the  variations  among  the  data 
sets.  The  results  of  this  process  are  illustrated  in  figure  7. 

compf  f  _in .  pro  -  The  purpose  of  this  routine  is  to  take  two  different  far-field  data  sets, 
display  both  of  them,  and  then  display  the  difference  between  them  in  a  similar  format.  As 
can  be  seen  in  figure  8,  the  difference  plot  is  rendered  using  a  specially  designed  color  table 
which  highlights  the  differences;  the  sign  of  the  difference  is  retained,  thereby  providing 
information  about  which  data  set  is  higher  than  the  other.  The  routine  also  allows  the  user 
to  select  or  experiment  with  the  dynamic  range  and  with  a  threshold  value  below  which 
differences  are  not  considered  important  and  can  be  blended  into  the  background.  This 
routine  calls  a  secondary  procedure  named  npf  f_m.pro;  the  PV-WAVE  code  for  this 
subordinate  procedure  is  not  included  but  is  functionally  equivalent  to  grabf  £ .  pro  and 
is  used  to  read  in  and  display  the  selected  data  sets. 
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£ig8  .pro  -  The  purpose  of  this  routine  is  to  use  wave- vector  processing  to  perform 
spatial  transforms  of  the  computed  normal  velocity  in  order  to  reveal  the  underlying  behavior 
of  elastic  waves  on  the  surface  of  the  scatterer.  It  produces  the  same  kind  of  helical  spectra 
that  have  been  used  in  holographic  work  with  experimental  data  to  associate  the  theoretical 
dispersion  curves  of  free  waves  on  an  infinite  shell  with  the  actual  motion  on  finite  closed 
bodies.  In  its  present  form,  it  is  only  appropriate  in  problems  for  which  there  is  a  cylindrical 
section  of  significant  length  compared  to  the  overall  dimensions  of  the  scatterer.  The  result 
of  this  routine  is  illustrated  in  figure  9,  in  which  the  bright  spots  indicate  waves  traveling  in 
one  axial  direction  with  a  particular  helical  wave  number.  Each  wave  number  is  an  (m,n) 
pair,  where  m  is  the  number  of  axial  half-wavelengths  along  the  length  of  the  cylinder  and 
n  is  the  number  of  the  wavelengths  that  fit  around  the  circumference  of  the  cylinder.  The 
pattern  of  these  spots  in  the  m,n  plane  has  the  appearance  of  a  “figure  8,”  which  gives  rise 
to  the  name  of  this  routine. 

window,  pro  -  This  is  a  sophisicated  use  of  the  PV-WAVE  software.  The  purpose  of 
this  routine  is  to  simulate  a  process  that  was  not  carried  out  in  the  computational  problem, 
but  which  controlled  the  accuracy  and  resolution  of  the  comparison  of  computed  results  with 
experimental  data.  Although  the  computation  in  the  Benchmark  Exercise  was  made  for  a 
free-field  steady-state  situation,  the  experimental  data  had  been  collected  in  a  small  tank 
requiring  a  short  time  window  to  be  applied.  This  routine  accepts  the  computational  data  in 
the  frequency  domain,  performs  a  series  of  fast  Fourier  transforms  to  achieve  a  real  time 
series  with  the  appropriate  sampling,  aligns  the  computed  time  series,  and  anplies  the 
windows  actually  used  in  the  experiment.  The  modified  time  series  is  then  transformed  back 
to  the  frequency  domain  for  comparison  with  the  experimental  data.  PV-WAVE  includes 
the  fast  Fc  rier  transform  and  other  signal  processing  functions  as  built-in  procedures/func¬ 
tions.  Figure  10  illustrates  the  end  result  of  using  window .  pro.  In  the  center  (figure  10b) 
are  experimental  results  in  the  normal  far-field  format,  showing  target  strength  in  color  as 
a  function  of  polar  angle  and  frequency.  On  the  right  (figure  10c)  the  computed  target 
strength  is  displayed  without  modification,  and  on  the  left  (figure  10a)  the  windowed 
computational  results  that  are  the  result  of  simulating  the  experiment  are  shown.  As  a  quality 
check ,  intermediate  time-domain  results  are  displayed  on  the  screen  while  the  routine  is 
running  through  the  data  by  angle,  as  shown  in  figure  11. 
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SUMMARY 


The  scattering  strength  of  elastic  objects  can  best  be  understood  by  displaying  its 
variations  in  creative  ways.  Futhermore,  visualization  of  such  large  multidimensional  data 
sets  provides  an  efficient  means  for  checking  as  well  as  understanding  or  analyzing  such  a 
function.  When  implemented  on  modem  graphical  workstations,  the  software  architecture 
and  routines  described  here  provide  a  powerful  environment  for  dealing  with  the  results  of 
computer-intensive  modeling  of  this  nature.  These  techniques  could  easily  be  modified  and 
employed  in  the  analysis  and  visualization  of  computational  or  'experimental  data  from  a 
wide  variety  of  physical  problems. 
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Figure  1.  Visualization  of  Benchmark  data  sets. 


Jobname  is  cla 

Comments:  CHIEF  with  2%  loss 

There  are  331  ka  values  from  0.200000  in  steps  cf  0.01 00000 
The  spectral  variable  is  ka 
There  are  1  ff  patterns 

24  areas  or  rings 
24  coordinate  values 
3  theta  incident  angles 

There  are  3  theta  inciderKe  angles  from  0.00000  in  steps  of  45.0000 

73  theta  observaUon  angles  in  ff  pattern  from  0.0  in  steps  of  2.50 

The  fluid  density  is  1 000.00 

The  sound  spe^  in  the  fluid  is  1500.00 

The  density  of  shell  is  7850.00 

The  Youngs  modulus  of  shell  is  2.00000e+11 

Poissors  ratio  for  the  shell  is  0.300000 

The  damping  factor  is  0.0200000 


Figure  2.  Example  output  from  rdhdr .  pro. 
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Max  FFP 
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Figure  3.  Complex  far-field  pressure. 
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Figure  4.  Normal  velocity,  clb,  90  incidence. 
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Figure  5.  Three-dimensional  representation  of  far-field  pressure  at  several  frequencies. 


Figure  6.  Axial  distribution  of  normal  velocity  for  various  frequencies. 


NPFF,  mag 


0:S0  0.75  130  1.85  2.40  2.95  330 

ka 

(b)  Envelope  of  selected  data 

Figure  7.  Data  set  variation  and  envelope  of  selected  data  sets. 
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Figure  8.  Comparison  and  differences  in  target  strength. 
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(a)  Windowed  (b)  Experiment  (c)  Computed  results 

Figure  10.  Comparison  of  experimental  results  with  windowed  and  unwindowed 

computational  results. 


0  500  1000  1500  2000  2500  0  500  1000  1500  2000  2500 


(a)  Reference  signal  (b)  Experimental  echo 


(c)  Computed  time  series  (d)  Windowed  time  series 

Figure  11.  Time  domain  processing  of  computational  data. 
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APPENDIX  A:  USTINGS  OF  PROGRAMS 


The  program  code  makes  use  of  a  number  of  built-in  PV-WAVE  functions  or  procedures. 
A  brief  functional  summary  of  the  PV-  WAVE  procedures  used  in  our  visualizations  is  given 
below.  These  summaries  are  extracted  with  permission  from  the  PV-WAVE  Technical 
Reference  Manual.  In  each  case,  the  syntax  of  the  Unction  or  procedure  is  given,  with  a 
brief  explanation  of  the  purpose  of  the  routine. 


Array  creation  routines 

BINDGEN(Dimi,...,Dimn) 

BYTARR(Dimi,...,Dimn) 

CDMPLEXARR(Dimi,...,Dimn) 

FlNDGEN(Dimi,...,Dimn) 

FLTARR(Dimi,...,Dimn) 

INDGEN(Dimi,...,Dimn) 

INTARR(Dimi,...,Dimn) 

LONARR(Dimi,...,Dimn) 

REPLICATE(Value,Dimi,...,Dmn) 

Array  manipulation  routines 

MAX(Array[,Max_Subscript]) 

MIN(Array[,Min_Subscript]) 

REBIN(Anray,Dimi,...,Dimn) 

REFORM(Array,Dimi,...,Dimn) 

SHlFl'(Array,Si,...,Sn) 

TRANSPOSE(Array) 

Data  conversion  routines 


Byte  array,  each  element  its  subscript. 
Returns  a  byte  vector  or  array. 

Complex  single-precision  vector  or  array. 
Floating  array,  each  element  its  subscript. 
Single-precision  floating  vector  or  array. 
Integer  array,  each  element  its  subscript. 
Returns  a  integer  vector  or  array. 

Returns  longword  integer  vector  or  array. 
Forms  new  array  filled  with  Value. 


Finds  the  maximum  element  of  an  array. 
Finds  the  minimum  element  of  an  array. 
Resamples  array  to  given  dimensions. 
Reformats  array  without  changing  contents. 
Shifts  the  elements  of  an  array. 

Transposes  an  array. 


BYTSCL(Array) 

FIX(Expr[,C)f£set[,Dimi,...,Dimn]]) 

STRING(Expri,...,Exprn) 


File  manipulation  routines 

CLOSE[,Uniti,...,Unitn] 
FREE_LUN,Uniti,...,Unitn) 
GET  LUN,Unit 
OPEMl,Unit,File 
OPENW,Unit,File 


Scales  and  converts  the  array  to  the  byte  type. 
Converts  parameter  to  integer  type. 

Converts  parameter  to  string  type. 


Qoses  one  or  more  files. 
Deallocates  one  or  more  files. 
Reserves  a  file  unit. 

Opens  a  file  for  reading  access  only. 
Opens  a  new  file  for  writing  access. 
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General  graphics  routines 

MOVIE,Images[,Rate] 

PLOTS, X[,Y[,Z]] 

SET_PLOTS,Device 
XYOUTS,X,Y, String 

General  Mathematical  functions 

ABS(X) 

CONJ(X) 

(Array, Direction) 

Image  display  routines 

TV,Image[,X,Y[,Channel]] 

TVSCL,Itnage[,X,Y[,Channel]] 

Image  processing  routines 

FFT(Array,Direction) 

Input  and  output  routines 

ASSOC(Unit,Array_structure[,Offset]) 
PRINT,Expri,...,Exprn  ) 
READ,Vari,...,Varn) 

READF,Unit,  Van,..,,  Yarn) 
READU,Unit,Vari,..,,Varn) 

Plotting  routines 

OPLOT,X[,Y] 

PLOT,X[,Y] 

Programming  routines 

HAK 

WAlT,Seconds 

String  processing  routines 

STRMID(Expr,First  Character,Length) 
STRTRIM(String,[flag]) 


Cycles  images  stored  in  3D  array. 
Draws  lines  (vectors)  and  points. 
Specifies  graphics  device. 

Sends  text  to  selected  graphics  device. 


Absolute  value. 
Complex  conjugate. 
Fast  Fourier  transform. 


Displays  an  image  on  the  display  screen. 
Scales  and  displays  an  image. 


Fast  Fourier  transform. 


Associates  variable  with  file  structure. 
Prints  to  the  standard  ouput  stream. 
Reads  from  the  standard  input  stream. 
Reads  from  the  specified  file  unit. 
Reads  unformatted  input. 


Plots  vector  arguments  over  old  axis. 
Plots  vector  arguments. 


"Hit  any  key  to  continue"  function. 
Delays  program  execution. 


Extracts  substring  of  string  expression. 
Removes  leading  and/or  trailing  blanks. 
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Transcendental  mathematical  functions 
AIX)G10(X) 

cospq 

SIN(X) 

Window  routines. 

WINDOW[,Window_Index]  Qeates  a  workstation  window. 

WSET[,Window_Index]  Selects  the  current  window. 


Base  10  logarithm. 
Cosine. 

Sine. 


A  list  of  the  Utility  and  Application  * .  pro  files  follows,  after  which  the  actual  programs 
and  example  output,  where  appropriate,  are  given.  In  the  programs,  any  line  or  portion  of  a 
line  beginning  with  a  semicolon  is  treated  as  a  comment. 

Utility  routines:  Commented  program  code  included 

convf fb.pro 
convbco.pro 
convsurb . pro 
rdhdr.pro 
grabff.pro 
grabsur.pro 
hardcopy, pro 
closelt.pro 

Application  routines:  Commented  program  code  and  example  output  included 

ff3d.pro 
plot sur. pro 
envel  pro 
compff_ni.pro 
figS.pro 
window. pro 
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;convffb.pro 

;La8t  modification:  20  April  92 


This  program  takes  a  binary  far-field  pressure  (.F£#)  file 
and  produces  an  equivalent  PV-Have  (.fw#)  file. 

After  this  conversion,  the  .ff#  is  not  needed  by  PV-Wave. 


10  ;  Select  the  data  set 

11 

12  fname=>8tring(replicate(32b,20) ) 

13  print, 'Input  job  name  -  no  extension' 

14  read, f name 

15  print, 'Input  incidence  angle' 

16  read, line 

17  case  iinc  of 

18  0: begin 

19  extension^' .ffO ' 

20  end 

21  45:begin 

22  extension^ ' . f f 4 ' 

23  end 

24  90:begin 

25  extension™' .ff9 ' 

26  end 

27  135: begin 

28  extension™' .ff3' 

29  end 

30  180:begin 

31  extension™' .ff8' 

32  end 

33  endcase 

34  fullname=fname-t-extension 

35 

36  get_lun,iunit 

37  openr , iunit , fullname , /f 77_unformatted 

3 8  jobname™string ( replicate ( 32b , 8 ) ) 

39  comment™string ( replicate ( 32b , 80 ) ) 

40  intrec™lonarr ( 1 0 ) 

41  rrecl*fltarr(10) 

42  rrec2™fltarr(10) 

43 

44  ;  Although  the  header  information  is  read  from  the  data  file, 

45  ;  it  is  not  stored.  The  header  information  from  the  .bco  file 

46  ;  is  used. 

47 

48  readu, iunit, jobname 

49  print ,' Reading  the  following  data  file:  ',fulln2une 

50  print, jobname 

51  readu , iunit , comment 

52  print, comment 

53  readu , iunit , intrec 

54  print, intrec 
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55  readu,iunit:,rrecl 

56  print, rrecl 

57  readu,iunit,rrec2 

58  print, rrec2 

59  print , intrec ( 0 ), '  ka  values  from  =  ', string (rrecl (0) ) 

60  print,'  in  steps  of  ', string (rrecl ( 1 ) ) 

61  print , intrec ( 1 ), '  different  patterns' 

62  print , intrec ( 2 ), '  theta  incidence  angles  from  ', rrecl (2) 

63  print,'  in  steps  of  ', rrecl (3) 

64 

65  ;  Next,  read  the  far-field  pressure  data 

66  ;  data  =  complexarr (nka,nharms,nobs,isymm) 

67  ;  assumes  coefficients  read  as  a  block 

68 

69  wector=coii^lexarr  ( intrec  ( 6 ) ) 

7  0  data=complexarr ( intrec ( 0 ) , intrec ( 5 ) , intrec ( 6 ) , int  rec ( 7 ) ) 

71  for  i=0, (intrec(O)-l)  do  begin 

72  for  1=0,  (intrec(5)-l )  do  begin 

73  for  m=0, (intrec(7)-l)  do  begin 

74  readu,iunit,wector 

75  data(i,l,*,m'=wec*-or 

76  endfor 

77  endfor 

78  endfor 

79  close, iunit 

80  print, 'Finioned  re»*ding  the  far-field  pressure' 

81 

82  ;  Next,  vrite  the  far-fiei.d  data  out  to  an  unformatted 

83  ;  PV-Wave  file 

84 

85  case  line  of 

86  C: begin 

8  7  extens ion= ' . f wO ' 

P8  end 

89  45:begin 

9  0  extens ion= ' . f w4 ' 

91  end 

92  90:begin 

9  3  extens ion= ' . f w9 ' 

94  end 

95  135:begin 

9  6  extens ion* ' . f w3 ' 

97  end 

98  180:begin 

9  9  extens ion* ' . f w8 ' 

100  end 

101  endcase 

102  outname=fname'fextension 

103 

104  get_lun, junit 

105  openw, junit, outname 

106  al=assoc(  junit, coii^lexarr(intrec(0)  ,intrec(S)  ,intrec(3)  ,$ 

107  intrec(7))) 

108  al(0)*data 
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109  £ree_lun,iunit 

110  £r«e_lun, junit 

111 

112  print Finished  writing  the  far-field  data' 

113  print, 'Finished  everything' 

114  end 
1 
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;convbco.pro 

jhtLBt  modification:  20  April  92 


This  program  takes  the  binary  .bco  fil&  and  produces  two 
ASCII  PV-Wave  files: 

.hdr  -  contains  all  the  parameters  of  the  job 
.wco  -  has  the  axial  and  radial  coordinates  of 
the  scatterer. 

After  this  conversion,  the  .bco  file  is  not 
needed  by  PV-Wave. 


14  ;  Select  the  data  set 

15 

16  f name^string (replicate (32b, 20 ) ) 

17  print, 'Input  job  name  -  no  extension' 

18  read, f name 

19 


20  ;  Open  the  file  with  header  and  coordinate  data  (.bco  file) 

21  ;  Then,  define  some  variables  and  read  the  header 

22 

23  get_lun , iunit 

24  openr, iunit, fname+' .bco' ,/f77_unformatted 

25  jobname=string(replicate(32b,8) ) 

26  comment=string( replicate (32b, 80) ) 

27  intrec>lonarr(10) 

28  rrecl«>fltarr(10) 

29  rrec2=fltarr(10) 

30  readu, iunit, jobname 

31  print ,' Reading  the  following  data  file:  ',fn2une 

32  print, jobname 

33  readu , iunit , comment 

34  print, comment 

35  readu , iunit , intrec 

36  print, intrec 

37  readu,iunit,rrecl 

38  print, rrecl 

39  readu, iunit, rrec2 

40  print,rrec2 

41  print , intrec ( 0 ), '  ka  values  from  =  ', string (rrecl (0) ) 

42  print,'  in  steps  of  ' ,string(rrecl ( 1 ) ) 

43  print , intrec ( 1 ), '  different  patterns' 

44  print,intrec(2) , '  theta  incidence  angles  from  ',rrecl(2) 

45  print,'  in  steps  of  ',rrecl(3) 

46 


47  ;  How  write  the  header  info  to  a  formatted  PV-Wave  file 

48 


49  get_lun, junit 

50  openw,  junit, fname-f' .hdr' 

51  printf, junit, jobname 

52  printf , junit, comment 

53  printf , junit , format- ' (10i5) ', intrec 

54  printf , junit, format-' (5gl0.4/5gl0.4) ' , rrecl 
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55  printf, junit,fonnat«' (5gl0.4/5gl0.4) ' ,rrec2 

56  free_lun, junit 

57  print,  'Finished  with  the  header  info' 

58 

59  ;  Next,  read  the  coordinates 

60 

6 1  coord^fltarr ( 2 , intrec ( 4 ) ) 

62  junk=f Itarr ( intrec ( 4 ) ) 

63  for  i>0,l  do  begin 

64  readu,iunit, junk 

65  coord(i,*)=junk 

66  endfor 

67  print, coord 

68  free_lun,iunit 

69 

70  ;  Now  write  the  coords  to  a  foarmatted  PV-Wave  file 

71 

72  get_lun, junit 

73  openw,  junit, fname-f' .wco' 

74  printf, junit, coord 

75  free_lun, junit 

76  print ,' Finished  with  the  coordinates' 

77 

78  print, 'Finished  everything' 

79  end 
1 
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;conv8urb.pro 

;La8t  modification:  20  April  92 


This  program  takes  a  binary  surface  data  file 

(either  surface  pressure,  .sp#,  or  normal  velocity,  .nv#) 

and  produces  an  equivalent,  binary  PV-Wave  file 

(either  .pw#  for  surface  pressure  or  .vw#  for  normal  velocity) 

After  this  conversion,  the  .nv#  or  .sp#  file  is  not  needed 

by  PV-Wave. 


13  ;  Select  the  data  type  and  set 

14 

15  fname«8tring(replicate(32b,20) ) 

16  surtype*'  ' 

17  print, 'Input  job  name  -  no  extension' 

18  read, f name 

19  print, 'Input  p  for  pressure  or  v  for  velocity' 

20  read,surtype 

21  case  surtype  of 

22  'p'tbegin 

23  stype^'.sp' 

24  otype='.pw' 

25  words^' surface  pressure' 

26  end 

27  'V' :begin 


30  words^ ' normal  velocity 

31  end 

32  endcase 

33  print, 'Input  incidence  angle' 

34  read, line 

35  case  line  of 

36  0: begin 

37  extension's  type-f '  0 ' 

3  8  oextens ion*otype+ ' 0 ' 

39  end 

40  45:begin 

4 1  extens  ion^s  type-l- '  4 ' 

4  2  oextens  ion>otype-f '  4 ' 

43  end 

44  90tbegin 

45  exten8ion>8type+'9' 

46  oextens  ion=otype-(- '  9 ' 

47  end 

48  135:begin 

49  extension^stype-f '  3 ' 

5  0  oextens  lon^otype-t- '  3 ' 

51  end 

52  180:begin 

53  extension»stype->-'8' 

5  4  oextens  ion*otype-i- '  8 ' 


55 

56 

57 

58 

59 

60 
61 
62 

63 

64 

65 

66 

67 

68 

69 

70 

71 

72 

73 

74 

75 

76 

77 

78 

79 

80 
81 
82 

83 

84 

85 

86 

87 

88 

89 

90 

91 

92 

93 

94 

95 

96 

97 

98 

99 
100 
101 
102 

103 

104 

105 

106 

107 

108 


end 

endcaee 

fullname^fname'fextension 
;  Open  the  file 

;  Then,  define  some  variables  and  read  the  header 
;  Header  is  read  from  the  data  file,  but  is  not 
;  written.  The  .hdr  file  is  generated  with 
;  convbco.pro 

get_lun , iunit 

openr , iunit , f ullname , /f 77_unformatted 

jobname>«string(replicate(32b,8) ) 

comiiient>string(  replicate  (32b,  80) ) 

intrec^lonarr (10) 

rreclsf Itarr ( 10 ) 

rrec2«f Itarr (10) 

readu , iunit , jobname 

print, 'Reading  the  following  data  file:  ',fullna]ne 

print , jobname 

readu , iunit , comment 

print , comment 

readu , iunit , intrec 

print , intrec 

readu , iunit , rrec 1 

print, rrecl 

readu , iunit , rrec2 

print, rrec2 

print , intrec ( 0 ) , '  ka  values  from  »  ' , string ( rrecl ( 0 ) ) 
print,'  in  steps  of  ', string (rrecl (1) ) 
print , intrec ( 1 ) , '  different  patterns ' 

print , intrec ( 2 ), '  theta  incidence  angles  from  ', rrecl (2) 
print,'  in  steps  of  ',rrecl(3) 

;  Next,  read  the  surface  data 

;  dumvector>complexarr ( nka , nharms , nareas , isymm) 

;  assumes  coefficients  read  as  a  block 

surf  >coiig>lexarr  ( intrec  ( 0 ) ,  intrec  ( 5 ) ,  intrec  ( 3 ) ,  intrec  ( 7 ) ) 
dumvector«complexarr ( intrec ( 3 ) ) 
for  i>0, (intrec(O)-l)  do  begin 

for  1«0, (intrec(5)-l)  do  begin 

for  m>0, (intrec(7)-l)  do  begin 
readu , iunit , dumvector 
surf ( i , 1 , * ,m) -dumvector 
endfor 
endfor 
endfor 

f ree_lun , iunit 

print , '  Finished  reading  the  '  -t-words 

;  Next,  write  the  surface  data  out  to  an  unformatted 
;  PV-Wave  file 
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109  outnuM»fnaine-(-oexten8ion 

110 

111  g«t_lun, junit 

112  openw, junltroutnane 

113  al>aasoc( junit,coi^lexarr(intrec(0) ,intrec(5) ,intrec(3) ,$ 

114  intrec(7))) 

115  al(0)«surf 

116  fraa_lun, junit 

117 

118  print, 'Finlahad  writing  the  '+%rardB 

119  print , ' Finished  everything ' 

120  end 
1 
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5 

6 

7 

8 
9 

10 

11 

12 

13 

14 

15 

16 

17 

18 

19 

20 
21 
22 

23 

24 

25 

26 

27 

28 

29 

30 

31 

32 

33 

34 

35 

36 

37 

38 

39 

40 

41 

42 

43 

44 

45 

46 

47 

48 

49 

50 

51 

52 

53 

54 


pro  rdhdr,filenajne>fna]ne 
;Laat  nodification:  20  April  92 


Thia  program  reads  header  from  a  data  set  with  the  same 
job  name 


comnon  header , jobname , comment , $ 

nka , npatts , nthetainc , nareaa , ncoords , nharms , nobs , isyimn, $ 
kastart , kainc , thetaincs tart , thetaincinc , phiincstart , $ 
phiincinc , thetaobsstart , thetaobsinc , $ 
rhof , cf , rhom, young , nu , eta 

;  First,  some  definitions 

jobname-string ( replicate ( 32b , 8 ) ) 
comment*atring( replicate (32b, 80) ) 
intrec^lonarr (10) 
rrecl>f Itarr ( 10 ) 
rrec2>°f  Itarr  (10) 


;  Now,  read  the  header  file 

get^lun, junit 

openr,  junlt,fname-f '  .hdr' 

readf, junit, jobname  ;  8  character  name  of  run 

readf, junit, comment  ;  80  characters  (20  words)  for  comments 

;  description,  etc. 

readf , junit , format* ' ( 1 0i5 ) ' , intrec 

readf , junit, format*' (5gl0.4/5gl0.4) ' ,rrecl 

readf , junit, format*' (5gl0.4/5gl0.4) ' ,rrec2 

f ree_lun , junit 

;  Give  the  header  info  some  more  meaningful  names 
;  First  the  integer  record 


nka*f ix ( intrec ( 0 ) ) 
npatts*f ix ( intrec ( 1 ) ) 
nthetainc*fix( intrec ( 2 ) ) 
nareas*fix ( intrec ( 3 ) ) 
ncoords*fix(intrec(4) ) 
nhanns*fix( intrec ( 5 ) ) 
nobs*f ix ( intrec ( 6 ) ) 
isymm*fix(intrec(7) ) 


;  number  of  ka  values 
;  number  of  ff  patterns 
;  number  of  theta  incident  angles 
;  number  of  areas  (rings) 

;  number  of  coordinates 
;  number  of  azimuthal  harmonics 
number  of  observation  angles  in  ff  patterns 
;  symmetry  flag,  1-yes,  2-no 


Now,  the  first  real  record 


kastart*rrecl ( 0 ) 
kainc*rrecl ( 1 ) 
thetaincstart*rrecl ( 2 ) 
thetaincinc*rrecl ( 3 ) 
phiincstart*rrecl ( 4 ) 
phlincinc*rrecl ( 5 ) 
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55  thataobsBtartsrrecl ( 6 ) 

5  6  thataobs inc»rrec 1(7) 

57 

58  ;  tha  rast  of  rracl  is  currently  not  used 

59  ;  Finally,  the  second  real  record 

60 


61  rhof ~rrec2 ( 0 ) 

62  cf>rrec2(l) 

63  rhoB^rrec2 ( 2 ) 

64  young*rrec2 ( 3 ) 

65  nu«rrec2(4) 

66  eta>rrec2(5) 


fluid  density 
sound  speed  in  fluid 
density  of  shell 
Young's  siodulus  of  shell 
Poisson's  ratio  of  shell 
damping  factor 


67 

68  ;  the  rest  of  rrec2  is  currently  not  used 

69  ;  Print  the  basic  information  about  this  data  set 

70 


71  print, 'Jobname  is  ',jobn2une 

72  print, 'Comments:  ', comment 

7  3  print , ' There  are  ' , nka , '  ka  values  from  ' , string ( kastart ) $ 

74  ,'  in  steps  of  ', string (kainc) 

75  print, 'The  spectral  variable  is  ka' 

76  print, 'There  are  ',npatts,'  ff  patterns' 

77  print,'  ',nareas,'  areas  or  rings' 

78  print,'  ',ncoords,'  coordinate  values' 

79  print,'  ' ,nthetainc, '  theta  incident  angles' 

80  print, 'There  are  ' , nthetainc , '  theta  incidence  angles  from  '$ 

81  ,thetaincstart, '  in  steps  of  ' ,thetaincinc 

82  print,'  ',nobs,$ 

83  '  theta  observation  angles  npff  pattern  from  ',$ 

84  thetaobsstart , '  in  steps  of  ' ,thetaobsinc 

85  print, 'The  fluid  density  is  ',rhof 

86  print, 'The  sound  speed  in  the  fluid  is  ',cf 

87  print, 'The  density  of  shell  is  ',rhom 

88  print, 'The  Youngs  modulus  of  shell  is  ', young 

89  print, 'Poissons  ratio  for  the  shell  is  ',nu 

90  print, 'The  dancing  factor  is  ',eta 

91 


92  end 


1 
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1  ;grabf£.pro 

2  ;Last  BKxlification:  20  ^rll  92 

3 

4  ;  — . — ____ — - - — - - - - - 

5  ;Thi8  progran  reads  in  PV>Wave  binary  far-field  file 

6  ;(.fw#)  placing  all  far->field  data  in  memory  for  use  as 

7  ;data. 

8  ; - 

9 

10  common  header, jobname, comment, $ 

1 1  nka , npatts , nthetainc , nareas , ncoords , nharms , nobs , isymm, kastart , $ 

12  kainc , thetaincstart , thetaincinc , phiincstart , phiincinc , $ 

1 3  thetaobsstart , the taobs inc , rhof , cf , rhom, young , nu , eta 

14 

15  jobname*'  ' 

16  pathname^'  ' 

17  print, 'input  jobname  (no  extension)  ' 

18  read, jobname 

19  print, 'input  pathname  including  trailing  /' 

20  read, pathname 

21  fullnamespathname-l- jobname 

22 

23  ;  read  vital  parameters  from  .hdr  file 

24 

25  rdhdr,filename>Bfullname 

26  print, 'Input  incidence  angle' 

27  read,ilnc 

28  case  iinc  of 

29  0:begin 

30  extension^' . fwO ' 

31  end 

32  45:begin 

33  extension^ ' . fw4 ' 

34  end 

35  90ibegin 

36  extension* ' . fw9 ' 

37  end 

38  135tbegin 

39  extension* ' . fw3 ' 

40  end 

41  180 (begin 

42  extension* ' . fw8 ' 

43  end 

44  else (print, 'Invalid  angle' 

45  endcase 

46  fullname*fullnametextension 

47 

48  ;  Open  the  file,  then  read  in  far-field  data 

49 

50  get_lun,iunit 

51  openr,lunit,fullname 

52  aa*assoc(iunit,coBqplexarr( nka, nharms, nobs , isymm) ) 

53  data*aa(0) 

54 
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55  ;  This  section  of  the  code  recombines  the  harmonic 

56  ;  con^nents  for  each  ka.  A  default  increment  in 

57  ;  phi  of  30  degrees  is  used. 

58  ;  Recombined  data  are  held  in  the  variable  npff. 

59 

60  deltaphi-BO 

61  nufflphi»l-(-fix(360./deltaphi) 

62  angles»deltaphi*( lpi/180. )*findgen(numphi) 

63  ivect»indgen(nharms) 

6  4  cosangle^transpose ( cos ( ivect#angles ) ) 

65  sinangle^transpose ( sin ( ivect#angles ) ) 

66  npff>complexarr(numphi,nobs,nka) 

67  for  kabin*0,nka-l  do  begin 

68  cteii^>BCOit^lexarr(numphi,nobs) 

69  if  isymm  eq  1  then  cteR^==coaangle#reform(data(kabin,*,*,0) ) 

70  if  isymm  eq  2  then  cteRq>scosangle#reform(data(kabin,*,*,0) )  $ 

71  +sinangle#reform(data(k£Lbin,*,*,l) ) 

72  npff  (*,*,kabin)=cteit^ 

7  3  endf or 

74  loadct,5 

75  f ree_lun , iunit 

76 

77  ;  This  section  of  code  produces  three  color  images 

78  ;  of  the  npff  data;  target  strength  in  db,  magnitude, 

79  ;  and  phase.  The  phi  observation  angle  is  fixed  at 

80  ;  0  degrees.  The  range  for  target  strength  is  [-30,30]. 

81  ;  The  range  for  magnitude  is  [0,max(mag) ] .  The  range 

82  ;  for  phase  is  [-180,180]. 

83  ;  In  each  case,  the  entire  theta  range  has  been 

84  ;  reconstructed  [0,360].  The  x-axis  corresponds  to 

85  ;  the  theta  observation  angle;  the  y-axis  to  the 

86  ;  ka  range. 

87 

88  xsz>2*nobs-l 

89  ysz>nka 

90  halfwayx:(numphi-l)/2 

91 

92  ;  Recombine  the  full  theta  observation  range 

93 

94  var»fltarr(xsz,ysz) 

95  var(0:nobs-l,*)=ab8(npff (0,*,*) ) 

96  var ( nobs : xsz-1, *) »abs ( reverse ( ref orm( npff ( halfway, l:nobs-l,*) ) ,1) ) 

97 

98  ;  Calculate  the  target  strength  in  dB 


This  section  of  code  produces  three  color  images 
of  the  npff  data;  target  strength  in  db,  magnitude, 
and  phase.  The  phi  observation  angle  is  fixed  at 
0  degrees.  The  range  for  target  strength  is  [-30,30]. 
The  range  for  magnitude  is  [0,max(mag) ] .  The  range 
for  phase  is  [-180,180]. 

In  each  case,  the  entire  theta  range  has  been 
reconstructed  [0,360].  The  x-axis  corresponds  to 
the  theta  observation  angle;  the  y-axis  to  the 
ka  range. 


100  ts*20*alogl0(var) 

101  tsimg»bytscl(ts,mins-30,max='30) 

102 

103  ;  Calculate  the  magnitude  and  phase 


105  magimg»byt8cl(var,min>0) 

106  pha8e>fltarr(xsz,ysz) 

107  phase( 0 tnobs-l, *)«( 180. / lpi)*atan (imaginary (npff(0,*,*) ) ,$ 

108  (float(npff(0, *,*))) ) 


1 09  phutt  ( nobs :  xsz>l ,  * )  ■‘reverse  ( ref  om(  ( 1 80 .  / 1  pi )  *$ 

110  atan(iinaglnary(np£f (halfway, l:nobB-l,*) ) ,£loat(npff (halfway, $ 

111  l:nobs-l, *)))),!) 

112  phaselng^bytscl (phase, inin«-180,max>180) 

113 

114  ;  Set  up  the  plotting  window 

115 

116  window, /free , xpoBs200 ,ypos«200 ,XBi2e‘>700 ,ysize=7 00 

117  tv,tsimg,100,100 

118  tv,Biagimg,270,100 

119  tv,phaseimg,440,100 

120  xyouts , / device ,100,450,' Target  Strength ' 

121  xyouts, /device, 270, 450, 'FF  Magnitude' 

122  xyouts, /device, 440, 450, 'FF  Phase' 

123  title* ' FAR-FIELO  DATA:  '+  comment 

124  xyouts, /device, 10 0,520, title 

125  title2*' Incidence  angle'  -t-  string(iinc) 

126  xyouts, /device, 150, 500, title2 

127  xyouts, /device, 400, 650, fullname 

128 

129  end 
1 
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1  ;grabaur.pro 

2  jLast  modification:  20  ^ril  92 

3 

4  , - 

5  ;Thi8  program  reads  in  a  PV-Wave  binary  normal  velocity  ( .vw#) 

6  ;or  surface  pressure ( .pw#)  file,  placing  all  surface  data 

7  ;in  memory  for  use  as  sur(nka,numphi,ncoords ) . 

8  ; - 

9 

10  common  header, jobname, comment, $ 

1 1  nka , npatts , nthetainc , nareas , ncoords , nharms , nobs , isymm , $ 

12  kastart , kainc , thetaincstart , thetaincinc , phiincs tart , $ 

13  phiincinc , thetaobsstart , thetaobsinc , rhof ,cf , rhom, young , $ 

14  nu,eta 

15 

16  fname*'  ' 

17  pathname='  ' 

18  print, 'input  jobname  (no  extension)' 

19  read, f name 

20  print, 'input  pathname  including  trailing  /' 

21  read, pathname 

22  fname>pathname+fname 

23 

24  ;  Read  essential  parameters  from  .hdr 

25 

26  rdhdr,filename=fname 

27  print, 'input  incidence  angle' 

28  read, line 

29  surtype*'  ' 

30  print, 'Do  you  want  surface  pressure(p)  or  velocity(v)?' 

31  read,surtype 

32  if  surtype  eg  'v'  then  begin 

33  title® 'NORMAL  VELOCITY  MAGNITUDE' 

34  8type='.vw' 

35  endif 

36  if  surtype  eg  'p'  then  begin 

37  title® 'SURFACE  PRESSURE  MAGNITUDE' 

38  8type*'.pw' 

39  endif 

40  case  iinc  of 

41  0: begin 

42  exten8ion®stype'f'0' 

43  end 

44  45:begin 

45  extension®stype+ ' 4 ' 

46  end 

47  90:begin 

4  8  extension®8type-f '  9 ' 

49  end 

50  135:begin 

51  extension®8type4- '  3 ' 

52  end 

53  180:begin 

5  4  extens ion®8 type+ ' 8 ' 
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55  end 

56  endcaae 

57  fullneune^fname-l-extension 

58 

59  ;  Open  the  file 

60 

61  get_lun , iunit 

62  openr ,  iunit ,  fullname ,  /£77_unf ormatted 

63  aa>a880C( iunit, coii^lexarr(nka,nhann8,ncoords,isyinm) ) 

64  data»aa(0) 

65  free_lun, iunit 

66 

67  ;  Thin  aection  of  the  code  recombinea  harmonica  for  each 

68  ;  ka  at  a  default  increment  in  phi  of  30  degreea. 

69 

70  deltaphi~30 

71  nus^hi>l+fix(360./deltaphi) 

72  anglea^deltaphi* ( lpi/180. )*findgen(numphi) 

73  ivect^indgen ( nharma ) 

7  4  coaangle-tranapoae ( coa ( ivect#angle8 ) ) 

75  8inangle=tran8poae(ain(ivect#angle8) ) 

7  6  8ur>*complexarr  ( nka ,  nun^hi ,  ncoorda ) 

77  ctempxcon^lexarr ( numphi , ncoorda ) 

73  for  kabin=0,nka-l  do  begin 

79  if  iaymm  eg  1  then  cteit^=co8angle#reform(data(kabin,*,*,0) ) 

80  if  iaymm  eg  2  then  cteiig)=co8angle#reform(data(kabin,*,*,0) )  $ 

81  +8inangle#reform(data(kabin,*,*,l) ) 

82  8ur(kabin,*,*)»cteii^ 

83  endfor 

84 

85  ;  Draw  labeled  picturea. 

86  ;  The  'Peraian  rug'  plot  ia  for  a  fixed  ka  (default  nka/2). 

87  ;  The  'atanding  wave'  plot  ia  for  a  fixed  phi  oba  angle  of  0  deg. 

88 

89  binno=nka/2 

90  phicut=0 

91  x»findgen(  nka  )*kainc-i-ka3  tart 

92 

93  ;  Draw  labeled  picture 

94 

95  print,min(ab8(8ur) ) ,max(ab8 (aur) ) 

96  print, 'input  maximum  value' 

97  read,maxval 

98  window, /free,xpo8=’100,ypo8=100,X8ize=900,yaize=650 

99  loadct,5 

100  y_8calesbindgen( 1 ,202 ) 

101  tv,congrid(y_8cale,30,331) ,800,70 

102  x4>[0,0] 

103  y»[0,maxval] 

104  plot, x4,y,X8tyle*4,y8tyle=l,po8=[ 800, 7 0,830, 401] ,/dev,$ 

105  /nodata, /noera8e,ytick3-7 

106  tv,byt8cl(congrid(ab8(reform(8ur(binno,*,*) ) ) ,20*nun^hi,288,$ 

107  /interp) ,min~0,max=maxval) ,90,70 

108  y-[0,l.] 
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109  x2-[0,360] 

110  plot , x2 , y , xsty le^l , ysty le»l , xtick8=6 , yt icks=4 , /nodata , $ 

111  /noerase ,  /dev ,  po8»  [90,70, 2  O^nusphi'f 9  0,358] 

112  xyout8 , /device ,210, 500, title, 8ize=2 

113  xyoute, /device, 300, 480,  'minval='-<-8tring(niin(ab8(8ur) ) ) 

114  xyoute, /device, 300, 465,  'niaxval='-t-8tring(max(ab8 (sur) ) ) 

115  tv,byt8cl(tran8po8e(congrid(re£onii(ab8(8ur(*,phicut,*) ) )  ,$ 

116  nka,288,/interp) ) ,mins0,max=maxval) ,425,70 

117  y-[0,l.] 

118  x3«[iiu.n(x)  ,inax(x)  ] 

119  plot , y , x3 , X8 ty le^l , ye ty le«l , xtick8=4 , yt ick8=l 1 , / dev , $ 

120  /nodata , /noeraee , po8=  [425,70,713, nka-l-7 0 ] 

121  xyoute, /device, 240, 570, 'incidence  angle  «  '+8tring(iinc) 

122  xyoute , /device ,240,620, £ul Inane , 8ize=2 

123  xyoute, /device, 120, 425, 'For  ka  =  '+Btrtrim(atring(x(binno) ) ,2) 

124  xyoute , /device ,450,425,' For  phibin  =  ' +8 trtrim ( e tr ing ( phicut ) , 2 ) 

125  xyoute, /device, 15 0,2 5, 'Phi' 

126  xyoute, /device, 4 5, 12 5, 'Normalized  axial  dietance' ,orientation=90 

127  xyoute, /device, 445,25, 'Normalized  axial  dietance' 

128  xyoute, /device, 3 80, 2 00, 'ka' ,orientation=90 

129 

130  end 
1 
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1  pro  hardcopy, FILENAME»fllename,BITS>bits,XSIZE=xsize,YSIZE=ysize 

2  ;La8t  modification:  24  Mar  92 

3 


1  pro  closeit 

2  ;La8t  modification:  21  April  92 

3 

4  - - 

5  ;ThiB  program  closes  the  hardcopy  file  and  sets  the  display 

6  ;bac]c  to  'x'  window,  replace  'x'  with  default  driver  on  system 

7  . - - 

8 

9  device, /close_file  ;  close  the  file,  so  it  can  be  sent  to 

10  ;  the  printer 

1 1  8et_plot , ' X '  ;  redirect  the  output  to  the  ' x '  window  driver 

12 

13  end 
1 
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1  ;ff3d.pro 

2  ;La(c  modification:  24  Mar  92 

3 


This  program  reads  in  PV-Wave  binary  far-field  file 
(.fw#),  generates  a  3-d  display  of  far-field  pressure  for 
user-defined  increments  in  ka,  and  displays  the  result  as  a 
8  ;  movie 


9  ; - 

10 

11  common  header, jobname, comment, $ 

12  nka , npatts , nthetainc , nareas , ncoords , nhanns , nobs , isymm, $ 

1 3  kastart , kainc , thetaincstzurt , thetaincinc , phiincstart , $ 

1 4  phiincinc , thetaobss tart , thetaobsinc , rhof , cf , rhom, young , nu , $ 

15  eta 

16 

17 

18  jobname>°'  ' 

19  pathnames'  ’ 

20  print, 'input  jobname  (no  extension)' 

21  read, jobname 

22  print, 'input  pathname  including  trailing  /' 

23  read, pathname 

24  fullnamespathnametjobname 

25  rdhdr,filenamesfullname 

26  print, 'Input  incidence  ^ngle' 

27  read,iinc 

28  case  iinc  of 

29  0:begin 

30  extensions ' . fwO ' 

31  end 

32  45:begin 

33  extensions ' . fw4 ' 

34  end 

35  90:begin 

36  extensions ' . fw9 ' 

37  end 

38  135:begin 

39  extensions ' . fw3 ' 

40  end 

41  180:begin 

42  extensions ' . fw8 ' 

43  end 

44  endcase 

45  fullnamesfullname+extension 

46 

47  get_lun , iunit 

48  openr, iunit, fullname 

49  aa«assoc( iunit, complexarr( nka, nharms, nobs, isymm) ) 

50  datasaa(O) 

51  free_lun, iunit 

52 

53  ;  This  section  of  the  code  recombines  harmonics  for  a 

54  ;  default  phi  increment  of  30  degrees.  Note:  a  finer 


55  ;  increment  may  be  desired  to  eliminate  plotting 

56  ;  artifacts . 

57  ;  define  some  variables 

58 

59  deltaphi=30 

60  nuiiqphiBl-l-fix(360./deltaphi) 

6 1  angles^deltaphi* ( 1 pi/ 1 8  0 . ) *  f indgen ( nus^hi ) 

62  ivectBindgen(nharms) 

6  3  cosangle^transpose ( cos ( ivect#angle8 ) ) 

64  sinangle-transpose ( sin ( ivect#angles ) ) 

65  npff»complexarr(nuitqphi,nobs,nka) 

66 

67  ;  form  npff  for  each  phi  angle 

68 

69  for  ]cabin>=0,nka-l  do  begin 

70  cteit^>coii^lexarr(numphi,nob8) 

71  if  isymm  eg  1  then  cteitg>=cosangle#  $ 

72  reform(data(kabin,*,*,0) ) 

73  if  isymm  eg  2  then  ctei[g>=co8angle#  $ 

74  reform(data(kabin,*,*,0) )  -l-sinangle#  $ 

75  ref orm( data ( kabin,*,*,l) ) 

76  npff (*,*,kabin)*cteng) 

77  endfor 

78 

79  ;  generate  frames  (default  of  5  ) 

80 

81  nframes=5 

82  8tep«fix(nka/nframe8) 

83  xp«fltarr(nob8,nurophi) 

84  yp>fltarr(nob8,numphi) 

85  zp>fltarr(nob8,nu]i^hi) 

8  6  indtheta>lndgen ( nobs ) 

8  7  indphi^indgen ( numphi ) 

88  unity«replicate(l,nun^hi) 

89 

90  ;  the  selection  of  the  radius,  a,  is  arbitrary,  but 

91  ;  it  does  affect  the  size  of  the  display 

92 

93  a-5. 

94  tang«lpi*indtheta/(nobs-l)  ;  generate  a  vector  of  theta's 

95  pang«2*lpi*indphi/(numphi-l)  ;  same  for  phi's 

96  zp«a*cos(tang)#unity 

97  factor>«a*sin(tang) 

98  xpBfactor#coB(pang) 

99  yp~factor#ain(pang) 

100  8h»bytarr ( nobs , numphi )  ;  this  holds  the  (bytscaled)  target 

101  ;  strength  values 

1 02  newlmg>bytarr  (370,365,  nf rames-l-l ) 

103 

104  window, 0,xpo8-350,ypo8»300,X8ize>:370,ysize»365 

105  loadct,5 

106  kaval>f indgen(nka) ^kainc+kastart 

107  for  j>0,nf reuses  do  begin  ;  freguency  loop 

108  k>j*8tep 
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109  8h*«bytscl (20*alogl0( transpose (aba (npf£(*,*,k) ) ) ) r$ 

110  min>-40,inax*40) 

111 

112  ;  the  next  routine  tends  to  choke  on  a  number  of  degenerate 

113  ;  polygons,  which  may  cause  some  artifacts,  but  the 

114  ;  routine  does  not  bomb  and  the  results  are  useful  as  is 

115 

116  8hade_surf _irr , sh*zp , sh*xp , 8h*yp , shades=sh , $ 

117  xrangeK[-1000, 1000] ,yrange«[-1000, 1000] ,zrange=[-1000, 1000] 

118  xyout8,0.4,0.8, 'ka  >  ' ■t-strtrim( string ( kaval( k) ))  ,$ 

119  size’ll , /normal 

120  wait, 0.0001 

121  newimg(*,*, j)«tvrd(0,0,370,365) 

122  endfor 

123 

124  ;  redisplay  the  frames  as  a  movie 

125 

126  movie , newimg , order^O 

127 

128  ;  make  a  hard  copy  of  any  frame  of  interest 

129 

130  ;  hardcopy ,  f  ile» '  clb9  3d#l  8 .  eps ' ,  xsize*3 .70,  ysize-^S .  6  5 

131  ;tv,newimg(*,*,18) 

132  ;closeit 

133 

134  end 
1 
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8 
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10 

11 


;plotsur.pro 

;Last  modification:  21  April  92 


This  program  makes  an  animated  multi-color  plot  showing  how 
surface  pressure  or  normal  velocity  vary  with  axial  node 
number  and  frequency.  It  assumes  that  grabsur.pro  has  been 
run  and  sur(nka,nuitq>hi,ncoords)  is  in  memory 


12  common  header, jobname, comment, $ 

1 3  nka , npatts , nthetainc , nareas , ncoords , nharms , nobs , isymm, $ 

14  kastart , kainc , thetaincstart , thetaincinc , phiincstart , phiincinc , $ 

1 5  thetaobsst art , that aobs inc , rhof , cf , rhom , young , nu , e t a 


16 


17  window, /free 

18  plot, sur(0, 0, *) ,ystyle=l,xstyle=l,yrange^[0,maxval] , /nodata,  $ 

19  xtitle«' Axial  bin  number' ,ytitle*' Magnitude' 

20  for  i*0,nka-l  do  begin 

21  oplot,abs(sur(i,0,*) ) ,color>i 

22  wait,. 02 

23  endfor 

24 

25 

26  end 
1 


A-25 


1 

2 

3 

4 

5 

6 

7 
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14 

15 

16 

17 

18 

19 
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22 

23 

24 


;«nv«'  ”»ro 

;Laat  loudification:  IS  May  92 

} - 

;Thls  program  generates  a  set  of  overlaid  spectral  plots  from 
; supposedly  equivalent  solutions  generated  by  different  codes, 

; allows  the  user  to  select  a  subset  of  the  solutions,  and 
; produces  a  mean  and  envelope  that  bound  the  variations  within 
;the  data  sets  selected. 

;Note:  This  program  assumes  the  data  to  be  co]iq>ared  has  been 
;read  in  prior  to  the  execution  of  envel.pro.  The  data  for 
;the  Benchmark  exercise  was  stored  in  variables  named  amagO 
; through  wmagO,  all  of  which  %rare  of  the  size  (nobs,nkas).  The  data 
;sets  were  reduced  by  the  selecting  the  phi  observation  angle  prior 
;to  execution. 

;The  example  given  was  produced  with  modified  versions  of  this 
; routine . 


; First  set  some  default  plotting  styles 

e 

§ 

lx. style*! 

ly. style«l 
lx.ticks>6 


25  ly.tieks*6 

26  ; 

27  id-'  ' 

28  print, 'Enter  index  of  theta  angle  desired' 

29  if  bmname  ne  '4a'  and  bmname  ne  '4b'  then  begin 

30  print, 'O-Oobs,  18«45obs,  36-90obs,  54-135obs,  72-180ob8' 

31  fact-2.5 

32  endif  else  begin 

33  print, '0-0obs,15-45obs,30-90ob8,45-135ob8,  60-180obs' 

34  fact-3.0 

35  endelse 


36  read, tangle 

37  ; 

38  print, 'Enter  number  of  data  sets  available' 

39  read,n8et8 

40  prlnt,nset8 

41  print, 'Enter  number  of  ka  values' 

42  read,nkas 

43  ; 

44  ;Set  up  x-axis  scale 

45  x-flndgen(nkaa)*.0l4.2 

46  if  (bmname  eq  '4a'  or  bmname  eq  '4b')  then  x-findgen(nkas) *.0094+ .397 

47  ; 

48  ;Detenalne  max(x)  as  legend  positioning  is  tied  to  this  value 

49  maxx-max(x) 

50  if (bmname  eq  'la')then  nka8-330  ;  special  case 

51  ; 

52  ;The  spectral  plots  are  saved  in  the  array  results. 

53  results-fltarr(nsets,nkas) 

54  results ( 0 , * ) -amagO ( tangle , 0 : nkas-1 ) 
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55  results ( 1 « * ) ■cmagO ( tangle , 0  snkas-l ) 

56  results ( 2 , * ) «snagO (tangle , 0  tnkas-1 ) 

57  results  ( 3 ,  * )  >ninnagO  ( tangle ,  0 :  nkas- 1 ) 

58  if  nsets  ge  5  then  result8(4,*)>niiiag0(tangle,0:nkas-l) 

59  if  nsets  ge  6  then  re8ult8(5,*)>:finag0(tangle,0:nka8-l) 

60  re8ult8ave»re8ult8 

61 

62  ;  The  y-axie  ie  autoacaled  to  l.l*max  of  the  data. 

63  maxdata-l . l^max ( reaulta ) 

64  ; 

65  ;  The  following  atatements  aet  up  line  typea  6 

66  ;  colora  and  the  legend  on  raw  data  plot. 

67 

68  Itype^intarr ( neeta ) 

69  ltype(0:3)>[0,0,0,4] 

7  0  colora^intarr ( naeta ) 

71  color8(0:3)>=[200,26,164,200] 

72  label>8trarr ( naeta ) 

73  label(0)*'l  -  axaar' 

74  label (!)■>' 2  -  chief' 

75  label(2)-'3  -  aara' 

76  label(3)«‘'4  -  waacat' 

77  if  naeta  ge  5  then  begin 

78  label(4)>‘'5  -  naahua' 

79  ltype(4)«4 

80  colora ( 4 )«26 

81  endif 

82  if  nsets  ge  6  then  begin 

83  ltype(5)-4 

84  colors ( 5 )-164 

85  label(5)-'6  -  fist' 

86  endif 

87  ; These  set  up  the  same  for  the  envelope  plot. 

88  lbl*[ 'average' , 'min  /  max'] 

89  ltyp-0 

90  clr-[164,200] 

91  ; 

92  window, /free,xaiz>625,y8iz«600,xpo8»0,ypo8-250 

93  loadct,12 

94  lx.title-'ka' 

95  Jy. title- 'NPFF,  mag' 

96  plot,x,re8Ults(0,*) ,yrange-[0,maxdata] , color-colors (0) 

97  for  i-l,nset8-l  do  oplot,x,re8ults(i,*) ,color-colors(i) ,$ 

98  lines ty le-ltype ( i ) 

99  legend, label, colors, ltype,p8ym, .7*maxx, . 87*maxdata,maxdata/40 . 

100  print, 'Input  problem  identifier;  e.g.  BMIB  45  INC' 

101  read, id 

102  id-strc<»ipres8(atrupcaae(id-(-'  ,  '+string(fix(fact*tangle)  )+•'  obs  ')) 

103  xyouts, .6, .90, /normal, id, color-200, size-1. 2 

104  xyouts, .6, .8 7, /normal, 'RAW  DATA' , color-200, size-1 .2 

105  ans-'y' 

106  while  ana  eq  'y'  do  begin 

107  maxdiff-0. 

108  setid-'  ' 
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rasults-resul'tsave 


109 

110 

111  ;Uaer  picks  how  many  and  which  data  sets  to  include  in 

112  ;tha  envelope.  The  user  can  chose  to  reprocess  the  envelope 

113 

114  print, 'how  many  data  sets  do  you  want  to  include' 

115  read,nunin 

116  if  numin  eq  nsets  then  setid-'ALL  DATA  SETS' 

117  if  numin  ne  nsets  then  begin 

118  setB*intarr(ntimin) 

119  print, 'which  data  sets  do  you  want  to  include', $ 

120  '1  -  nsets  possible' 

121  read, seta 

122  print, 'you  have  asked  to  include  sets', sets 

123  for  i*0,nufflin-l  do  setid-aetid-fstring(sets(i) ) 

124  for  i«0,numin-l  do  seta(i)>sets(i)-l 

125  results^results ( sets , * ) 

126  endif 

127  average»avg(results,0) 

128  window, /free, xsiz=625,ysiz«600,xpos«625,ypos»250 

129  mini’-min(  results) 

130  ]naxdum»fltarr(nkas) 

131  mindum>fltarr(nkas) 

132  for  i»0,nkas-l  do  begin 

133  maxdum( i)=Biax (results (*, i) ) 

134  mindum(  i)>inln(  results  (*,  i) ) 

135  endfor 

136  Biaxd«1.2*inax(aiaxdum) 

137  plot,x,  (average) ,  /nodata, yrange>[0,iiiaxdat a] 

138  oplot,x, (average) , color-164 

139  oplot , X ,  (maxdum) 

140  oplot, X,  (mindum) 

141  legend, lbl,clr,ltyp,psym, .7*inaxx, .8*maxdata, .3 

1 42  xyouts , . 6 , . 85 , /normal , BtrcompreBB ( setid ) , color-2  00, siz-1 . 2 

143  xyouts, .6, .90, /normal, id, color-2 00, siz-l .2 

144  maxdiff-iDax(iaaxdum-mindum) 

145  print, 'Maximum  difference  -  ',inaxdiff,'  at  ',$ 

146  x(where(maxdum-mindum  eq  maxdiff ) ) 

147  ; oplot, x,n tru th ( 0, *) ,color-100 

148  print, 'Do  you  want  to  reprocess,  y  or  n' 

149  read,ans 

150  endwhile 

151  end 


1 

2 

3 

4 

5 

6 

7 

8 
9 

10 

11 

12 

13 

14 

15 

16 

17 

18 

19 

20 
21 
22 

23 

24 

25 

26 

27 

28 

29 

30 

31 

32 

33 

34 

35 

36 

37 

38 

39 

40 

41 

42 

43 

44 

45 

46 

47 

48 

49 

50 

51 

52 

53 

54 


pro  coii^ff_m 

9 

;LaBt  nodificatlon:  19  May  1992 


This  procedure  plots  the  target  strength  image  for  two  codes  and 
plots  the  difference  between  the  images. 

The  difference  is  taken  between  dB  values  which  have  first  been 
subjected  to  a  user-specified  threshold. 

The  user  is  allowed  interactive  control  over  the  dynamic  range 
of  the  images  and  the  difference  image. 

Note:  The  procedure  assumes  npff_m.pro  has  been  con^iled. 


common  header,  jobname,  comment,  $ 

nka , npatts , nthetainc , nareas , ncoords , nharms , nobs , isymm , $ 
kastart , kainc , thetaincstart , thetaincinc , phiincstart , phiincinc , $ 
thetaobsstart , thetaobs inc , rhof , cf , rhom, young , nu , eta 

e 

9 

common  sharel,  npff,  nun^hi,  diff 

9 

;  Get  the  first  set  of  data 

e 

9 

print , ' Prompting  for  the  threshold  value  in  db ' 
read , thresh_db 

print ,' Prompting  for  dynamic  range  desired  for  difference' 
read , dmin , dmax 

print ,' Prompting  for  the  first  far-field  data  file.' 
npff_m 

npff l*abs( npff ) 

;npffl(where(abs(npffl)  le  thresh) )«thresh 
nkal>nka 

e 

9 

print ,' Proo^ting  for  the  second  far-field  data  file.' 
npff_m 

npf  f  2»2d>s  ( npf  f ) 

; npff2 (where (abs(npff2)  le  thresh) ) -thresh 

• 

9 

xsz*2*nob8-l 

ysz>nka 

9 

yszl>nkal 

; Calculate  the  target  strengths 

e 

9 

ansd« ' n ' 
anst> ' n ' 
ans> ' y ' 
repeat  begin 

print , ' Input  the  index  of  the  phi  angle  (0  to  6 ) . ' 
read,philnc 

halfway»(nus^hi-l)/2  +  phiinc 

varl^f Itarr ( xsz ,yszl ) 

varl ( 0 :nob8-l , * )»(npf f 1 (phiinc, * , * ) ) 
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var 1 ( nobs : xsz-l , * ) = ( reverse ( ref onn( npf f 1 ( halfway , $ 

1 :nobs-l,*) ) , 1) ) 

tsl«20*alogl0(varl) 

var2»f Itarr ( xs z , ys z ) 

var2(0:nobs-l,*)»(npff2(phiinc,*,*) ) 

var2 ( nobs : xsz-1 , * ) - ( reverse ( ref onn( npf f 2 ( halfway , $ 

l:nobs-l,*) ) ,1) ) 

ts2»20*alogl0 ( var2 ) 

loadct,5 

print, 'max  values' ,max(tsl) ,max(ts2) , '  min  values ',$ 
min(tsl) ,min(ts2) 

print , ' input  dynamic  range  for  ts  pictures ' 

read , mints , maxt s 

top-202 

print, 'input  number  of  colors  on  system' 
read,topval 

window, 0 , xsize-2*xsz ,ysize-2*ysz , xpos-1 00 ,ypos=200 , $ 
title- 'Set  1' 

tv , rebin ( bytscl ( ts 1 ,min=mints ,max-maxts , top-topval ) , $ 
2*xsz,2*yszl) 

window , 1 , xs i ze=2 *xs z , ys ize=2 *ys z , xpos= 8 0 0 , ypoa=2 0 0 , $ 
title- 'Set  2' 

tv , rebin ( bytscl ( ts2 ,min-mints ,max-maxts , top-topval ) , $ 
2*xsz,2*ysz) 

window, 5, X8ize-100,ysiz-300,xpos-1125,ypos-550 
y_8cale-bindgen ( 1 , topval ) 
tv,congrid(y_8cale,30,250) ,45,30 
x-[0,0] 

y- [  mint  8 ,  maxts  ] 

plot ,x,y ,X8tyle-4 ,y8tyle=l ,pos- [45,30,75,280], /device, $ 

/nodata , /noerase , yticks-4 
print , ' hit  when  ready ' 
hak 

repeat  begin 

tsl (where (t8l  le  thresh_db) )-thresh_db 
ts2 (where ( ts2  le  thresh_db) ) -thresh_db 
diff-t8l-ta2 
repeat  begin 
loadct , 1 1 
;8tretch, 0,255 

window, 2 ,x8ize-2*xsz,y8ize-2*ysz ,xpos-450 ,ypos-200 , $ 
title- ' Difference ' 

tv , rebin ( by tael ( dif f ,min-dmin , max-dmax , top-topval ) , $ 
2*xsz,2*y8zl) 

window, 6 , xsiz-lOO ,ysiz-300 ,xpo8-1125 ,ypos-200 

tv,congrid(y_scale,30,250) ,45,30 

x-[0,0] 

y-[dmin,dmax] 

plot , X , y , xs ty le-4 , ysyle-l , po8- [45,30,75,280], /device , $ 
/nodata , /noerase , yticks-4 , color-101 
print, 'Do  you  want  a  new  dynamic  range  for  the  difference* 
read,ansd 

if(an8d  eg  'y')  then  begin 

print, 'Input  min  and  max  desired' 
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109  read,dmin,dmax 

110  endif 

111  endrep  until  ansd  eq  'n' 

112  print, 'Do  you  want  another  threshold  value?  (y/n) ' 

113  read, ana t 

114  if(anst  eq  'y')  then  begin 

115  print, 'Input  new  threshold  value' 

116  read,thresh_db 

117  endif 

118  endrep  until  anat  eq  'n' 

119  print, 'Do  you  want  another  phi  angle?  (y/n)' 

120  read, ana 

121  endrep  until  ans  eq  'n' 

122  end 
1 
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1  ;fig8.pro 

2  ;Laat  modification:  20  Jul  92 

3 

4  ; - 

5  ;This  program  generates  and  displays  color  surfaces 

6  ; which  portray  target  strength  as  a  function  of  frequency 

7  ; and  azimuth 

8  ; - 

9 

10  common  header, jobname, comment, $ 

1 1  nka , npatts , nthetainc , nareas , ncoords , nharms , nobs , isymm, $ 

12  kastart ,kainc, thetaincstart,thetaincinc,phiincstart , $ 

1 3  phiincinc , thetaobsstart , thetaobsinc , $ 

14  rhof ,cf,rhom, young, nu, eta 

15 

16  datapath= ' /data/schenck/dset2/ ' 

17  jobname='  ' 

18  print, 'input  runname  (without  extension' 

19  print,'  assumed  to  be  in  /data/schenck/d8et2 ) ' 

20  read,jobname 

21 

22  fullname=datapath-t-jobname 

23  rdhdr,filenamesfullname 

24 

25  print, 'Input  incidence  angle' 

26  read, line 

27  case  iinc  of 

28  0:begin 

29  extensions^' .vwO' 

30  end 

31  45:begin 

32  extensions' .vw4 ' 

33  end 

34  90:begin 

35  extensions' .vw9 ' 

36  end 

37  135:begin 

38  extensions' .vw3 ' 

39  end 

40  180:begin 

4 1  extens ion* ' . vw8 ' 

42  end 

43  endcase 

44  fullnamesfullname+extension 

45 

46  get_lun , iunit 

47  openr, iunit, fullname 

4 8  aasassoc ( iunit , complexarr ( nka , nharms , ncoords , isymm) ) 

49  velsaa(O) 

50  free_lun , iunit 

51 

52  window, 2, xsizes50,ysizes512,xpossi50,yposs300 

53  loadctfll 

54  cbarsbytarr(50,512) 


55  for  i*0,255  do  cbar(*,i*2:i*2+l)=i 

56  tVfCbar 

57 

58  nendnodes=12  ;  presently  not  in  header  file 

59  npt8«narea8-2*nendnode8 

60  infact»5  ;  multiplication  factor  for  display 

61  nfact<*5  ;  multiplication  factor  for  display 

62  window, 4 ,xpos=250 ,ypos=375 ,xsize»mf act^npts , $ 
S3  y8ize=nfact*(2*nharma-l) 

64 

65  char^'g'  ;  a  junk  character 

66 


67  for  1=0,0  do  begin 

68 

f=9 0+10*1 

if  isymm  eg  2  then  begin 

vf =ref orm( vel ( f , * , nendnodes : nendnodes+npts-l , * ) ) 
vh=complexarr ( 2*nharms-l , npts ) 
for  h=-(nharms-l) ,nharms-l  do  begin 
if  h  eg  0  then  begin 

vh(h+nharms-l,*)=reform(vf (h,*,0) ) 
endif 

if  h  gt  0  then  begin 

vh ( h+nharms-1 . * ) =0 . 5*ref orm( vf ( h , * , 0 ) ) -com- 


69 

70 

71 

72 

73 

74 

75 

76 

77 

78 

plex(0,l)*$ 

79 

80 
81 
82 

plex(0,l)*$ 

83 


reform(vf (h,*,l) ) 
endif 

if  h  It  0  then  begin 

vh ( h+nharms-1 , * ) =0 . 5*ref orm( vf ( -h, * , 0 ) ) -com- 


reform(vf (-h,*,l) ) 

84  endif 

85  endfor 

86  endif  else  begin  ;  to  allow  cosine  only  representation 

87  vf=reform(vel(f ,*,nendnode8:nendnodes+npts-l) ) 

88  vh=con^lexarr(2*nhanns-l,npt8) 

89  for  h=-(nharms-l) ,nhanns-l  do  begin 

90  if  h  eg  0  then  begin 

91  vh ( h+nharms-1,*) =ref orm ( vf (h,*) ) 

92  endif 

93  if  h  gt  0  then  begin 

94  vh(h+nharms-l,*)=0.5*reform(vf (h,*) ) 

95  endif 

96  if  h  It  0  then  begin 

97  vh(h+nharmfl-l,*)=0.5*reform(vf (-h,*) ) 

98  endif 

99  endfor 

100  endelse 

101  npts=47  ;  ###  to  force  npts  to  be  odd,  not  general! 

102  vh=vh(*,0:46)  ;  same 

103 

104  kmat=coRqplexarr(2*nharms-l,npt8) 

105  8ym=coiiq>lexarr(2*nharms-l,npts) 

106  middle- ( npts- 1 ) /2 
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107  for  i»0,2*(nhanns-l)  do  begin 

108  vectO=re£orm(vh(i,*) ) 

109  Jcmat(i,*)»fft(vectO,l) 

110  8ym(i, middle )>kmat(i,0) 

111  unitv^Hioii^lexarr  (middle-i-1 ) 

112  unitv(0)«coii^lex(l,0) 

113  unitv(l)>coii^lex(cos(  lpi*(npts-l)/npts)  ,$ 

114  sin( lpi*(npts-l)/npts) ) 

115  for  j>l, middle  do  begin 

116  unitv( j)=unitv(l)*unitv( j-1) 

117  i »middle+ j ) =unitv ( j ) *kmat ( i , npts- j ) 

118  8ym( i, middle- j)=conj(unitv(j) )*kmat(i, j) 

119  endfor 

120  endfor 

121  kmat>8ym 

122  new«transpo8e(kmat) 

123  8ign*=float  (new) /ab8(  float  (new) ) 

124  newer«8ign*ab8 (float (new) ) 

125  tv,rebin(bytscl(newer,min»-25,max»25,top«254) ,mfact*npts,$ 

126  nfact*(2*nharm8-l) ) 

127 

128  xyout8, 0.05, 0.05, 'freq»'+8trtrim( string ( (f+l)/2. ) ,2) ,$ 

129  /normal,colorB254 

130 

131  print, 'Frequency  bin  =  ',£ 

132  wait, 0.5 

133  char>get_kbrd(0) 

134  if  char  ne  ' '  then  begin 

135  print, 'Bit  any  key  to  continue' 

136  hak 

137  endif 

138 

139  endfor 

140  end 
1 
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1  ;  window. pro 

2  ;  Last  modification:  9  Jun  92 

3 

4  • - 

5  ;  This  program  takes  model  data  in  the  frequency  domain  and  performs 

6  ;  the  windowing  prescribed  by  URL,  then  presents  a  comparison  of  the 

7  ;  model  results  with  the  experimental  data  in  the  frequency  domain 

8  ; - 

9 

10  common  header, jobname, comment, $ 

1 1  nka , npatts , nthetainc , nareas ,ncoordB , nharms , nobs , isymm, $ 

12  kastart , kainc , thetaincstart , thetaincinc , phiincst^u:t , phiincinc , $ 

1 3  thetaobss tart , thetaobs inc , rhof , cf , rhom, young , nu , eta 

14 

15  ;  Select  the  data  set 

16 

17  ;  Change  pathname  data  file  n2unes  as  required 

18  pid*'  ' 

19  print, 'input  pid,  eg  4a  4b' 

20  read, pid 

21  pathname^ ' /data/schenck/dset2/ ' 

22  ;pathname‘>'  ' 

23  ; print , ' input  pathname ' 

24  ; read , pathname 

25  jobname='  ' 

26  print ,' Input  jobname  (no  extension)  ' 

27  read, jobname 

28  print, 'Input  incidence  angle  in  degrees' 

29  read,iinc 

30  case  pid  of 


31 

' 4a' :begin 

32 

refname*'bm4aref .daw'  ;reference  time  series 

33 

case  iinc  of 

34 

0 : begin 

35 

extension*' .fwO' 

;  con^utational  data 

36 

cdat* ' cleanaO . dat ' 

;  cleaning  window  data 

37 

shifty=-0  ;  to  account  for  orientation  of  target 

38 

;  in  NRL  tank 

39 

f name* ' bm4a0 . daw ' 

;  experimental  echo  time  series 

40 

end 

41 

42 

45:begin 

43 

extension* ' . f w4 ' 

44 

cdat* ' c leana4 . dat ' 

45 

shifty*- 15 

46 

f  name* ' bm4  a4  5 . daw ' 

47 

end 

48 

49 

50 

51 

52 

53 

54 


90:begin 

extension*' .fw9 ' 
cdat= ' cleana9 . dat ' 
shifty=-30 
f name* ' bm4 a9  0 . daw ' 
end 


endcase 
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80 
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86 


end 

'4b' : begin 

re£name> '  bm4bre£ .  daw' 
case  line  o£ 

0:begin 

extension^' .£w0' 
cdat» ' cleanbO . dat ' 
shi£ty»-0 
£naines '  bm4b0 .  daw ' 
end 

45:begin 

extension’s ' .  £w4 ' 
cdat= ' cleanb4 . dat ' 
ahi£ty*-15 
inastes: '  bm4b4  5 .  daw ' 
end 

90:begin 

extension- ' . £w9 ' 
edat’s '  cleanb9 .  dat ' 
ahi£ty=-30 
£naine= '  bm4b9  0 .  daw ' 
end 
endcase 

end 

else  t print , ' Incorrect  jobnane ' 
endcase 

iullnane^pathnaine-l- jobnane 


87  ;  Read  the  header  and  open  the  £ile 

88  rdhdr,£ilenane»£ullname 

89  £ullname=£ullnane-i-extension 


90  get_lun,iunit 

91  openr,iunit,£ullname 

92  aa>assoc(iunit,coii^lexarr(nka,nhams,nobs,isynm) ) 

93  data»aa(0) 

94  £ree_lun,iunit 

95  print, 'Finished  reading  nodel  data  set' 

96 


97  ;  this  next  part  could  be  made  simpler, 

98  ;  since  we  only  want  phicut^O 

99 


100  deltaphi>180 

101  nu3i^hi-l+£ix(360./deltaphi) 

102  angles>daltaphi*(  ipi/180.  )*£indgen(nu]i^hi) 

103  ivect*>indgen  ( nharms ) 

104  cosangle-transpoae(cos(ivect#angles) ) 

1 05  sinangle>transpose ( sin  < ivect#angles ) ) 

106  np££>coiiq^lexarr(nuiiqphi,nobs,nka) 

107  £or  kabin*0,nka-l  do  begin 

108  ctemp*complexarr ( nua^hi , nobs ) 
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109  if  isyimn  eq  1  then  cteiiq>«co8angle#refonn(data(kabin,*,*,0) ) 

110  If  iaymm  eq  2  then  cteii)pacosangle#refonn(data(k2d>in,*,*,0) )  $ 

111  -t-8inangle#refonD(data(k2d}ln,*,*,l) ) 

112  npff (*,*,kabin)=cten^ 

113  endfor 

114 

115  ;  Calculate  the  full  360  degree  pattern  of  computed  target  $ 

116  ;  strength 

117 

118  xsz*:2*nobs-l 

119  ysz^nka 

120  halfway<‘(nu]nphl-l)/2 

121  var»coiiq>lexarr(xsz,ysz) 

122  var(0;nob8-l,*)=npff (0,*,*)  ;  phicut  =  0 

123  var(nobs:xsz-l,*)=reverse(reform(npff (halfway, l:nobs-l,*) ) ,1) 

124  print, 'Full  azimuthal  pattern  constructed' 

125 

126  window , 5 , xs ize=9  0  0 , ys ize=4  00 , xpo8=2  00, tit le=$ 

127  '5. Target  strength  vs  angle  and  frequency' 

128  loadct,5 

129  bytc»bytscl(20*alogl0(aLbs(var) )  ,mins-30,iiuuc=30,$ 

130  top*ld.n_color8-l) 

131  tv, rebin(bytc, 2*121, 331) ,600,35 

132  xyouts, 0.72, 0.95, 'Unwindowed  model ', /normal 

133 

134  add«42  ;  adjust  lowest  non->zero  bin  of  computed  data 

135  ;  (fir8tbin*binwidth=1.7kHz) 

136 

137  ;  First,  we'll  work  with  the  reference  signal 

138 

139  ;  Read  in  ttie  reference  (incident)  time  series 

140 

141  refname^'  /data/schenck/expt4/'-i-refname 

142  get_lun,iunit 

143  openr , iunit, ref name, /f77_unformatted 

144  bb-a88oc(iunit,fltarr(4096) ) 

145  ref=bb(0)  ;reference  (incident)  time  series 

146  free_lun, iunit 

147  tape_info=ref (4080:4095) 

148  print, tape_info 

149  ref (4080:4095)=0.  ;  replace  last  16  bins  with  zero 

150 

151  ref>8hift(ref ,380) 

152  window, l,xsize=400,ysize>200,xpos=150,ypos=7 00,$ 

153  title*' 1 .Reference  signal' 

154  plot,ref (0:2500) ,yrange*[-25, 25] ,ystyle=l  ;  plot  the  time  series 

155  tlength*12423  ;  to  achieve  2  usee  san^le  period 

156  ;  in  the  coiiq>uted  time  series 

157  rpad*replicate(0. 0,tlength) 

158  rpad(0:4095)*ref 

159  cpad*fft(rpad,-l) 

160  magpinc*avg ( abs ( epad ( add : addtnka- 1 ) ) ) 

161 

162  ;;;  change  pathnames  if  required 


163  fullnane^' /data/Bchenck/expt4/ '+fname 

164  get_lun,iunit 

165  op«nr,lunit,fullnaina, /f 77_un£onnatted 

166  aa«assoc(iunit, fltarr( 121, 4096 ) ) 

167  acho»aa(0) 

168  £ree_lun,iunit 

169  print, 'Finished  reading  the  experimental  data' 

170 

171  ;  The  correction  £actor  accounts  £or  the  source, 

172  ;  target,  receiver  geometry 

173  geocorrect>2*2*sqrt (2. 18)/ (2.239*2.1 8*0. 0254*sqrt(l .96) ) 

174  ech>geocorrect*echo 

175 

176  window, 2, xsize>400,ysizea200,xpos«150,ypoB«470,$ 

177  title*' 2 .Experimental  echo' 

178  loadct,5 

179  window, 3, xsize*400,ysize*200,xpos*150,ypos*10,$ 

180  title*' 3. Computed  time  series' 

181  window, 4 ,xsize*400 ,y8ize*200,xpos*150 ,ypos*240 , $ 

182  title* '4. Windowed  time  series' 

183 

184  ££sur£*coiig>lexarr(121,331) 

185  modspec*coiiq>lexarr(  12 1,331) 

186 

187  ;  Mow  develop  cosine  squared  window  £or  each  observation  angle 

188  ;  this  ln£ormation  comes  £rom  Brian  Houston's  memo  re  £rom 

189  :  3600  to  3700  data 

190 

191  cdat*'-/prob4/'+cdat 

192  get_lun,kunit 

193  openr,kunit,cdat 

194  cleanarr*intarr(4,121) 

195  read£,kunlt,cleanarr 

196  £ree_lun,kunit 

197 

198  print, cleanarr 

199 

200  ;  Now  we'll  start  the  large  loop  over  obsexrvation  angle 

201 

202  £or  j*0,120  do  begin 

203 

204  k*-( j+Bhi£ty) 

205  i£  k  It  0  then  k«k+121 

206 

207  ;  Next,  let's  process  the  experimental  echo 

208 

209  print, echo(k, 4094)  ;  print  angles  of  observation 

210  ech(k, 4080: 4095 )*0.  ;  zero  the  last  sixteen  bins 

211 

212  longecho*replicate ( 0 . ,t length) 

213  longecho(0:4095)*re£orm(ech(k,*) ) 

214  back8pec*££t(longecho,-l)  ;  convert  to  £requency  domain 

215  ££8ur£ ( j , * ) *back8pec ( add : addtnka- 1 ) 

216 
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217  %r8et,2 

218  plot, echCk, 0:2500) ,yrange>[-10. ,10. ] ,ystyle>l 

219 

220  ;  Generate  and  overplot  the  cleaning  window 

221 

222  frontstart«cleanarr(0,k) 

223  frontend>cleanarr(l,k) 

224  back8tart>cleanarr(2,k) 

225  backend«cleanarr(3,k) 

226  flength^frontend-frontstarttl 

227  blength>backend-backstart-(-l 

228  piover2>Ipi/2. 

229  frontw«(8in(piover2*findgen(£length)/float(flength-l) ) )^2 

230  backw«( cos (piover2*£indgen(blength) /float (blength-1) ) )^2 

231  cleanw>replicate(0.0,4096) 

232  cleanw( f rontstart : f rontend) rontw 

233  cleanw(  frontend-l-1  tbackstut-l )  >1 . 0 

234  cleanw(backstart: backend )>backw 

235 

236  longclean>replicate ( 0 . , tlength) 

237  longclean(0:4095)-cleanw 

238  oplot,5.*longclean(0:2500) ,color~120 

239  wait, 0.02 

240 

241  ;  Next,  we  will  process  the  coaqputed  results 

242 

243  bsf£«re£om(var(  j,*)  )*cpad(add:add+nka-l)  ;  positive  freqs 

244  part*replicate( 0.0, tlength) 

245  padded>complex( part, part)  ;  a  complex  atrray  of  zeroes 

246  padded (0‘t-add:add't-nka-l)>bsff  ;  fill  in  the  negative  frequencies 

247  rev^onj  (reverse (bsff ) )  ;  ensures  a  real  time  series 

248 

249  ;  fill  in  the  negative  freqs 

250  padded(tlength-nkatl-add:tlength-add)>rev(0:nka-l)  ; 

251  tseries*f ft (padded,!)  ;  take  the  FFT 

252  aK}dtO>Bfloat(tseries) 

253 

254  W8et,3 

255  plot,modt0,yrange>[-10.0,10.0] ,ystyle>l 

256  oplot, imaginary (tseries) ,color*120  ;  this  should  be  zero 

257  wait, 0.02 

258 

259  ;  Next,  apply  the  cleaning  window  to  the  coiiq>uted  results 

260 

261  newc*longclean*modt0  ;  apply  the  window  in  time  domain 

262  newspec»fft(newc,-l )  ;  compute  the  windowed  model  spectrum 

263  modspec ( j , * ) "newspec ( add : addtnka-l ) 

264 

265  wset,4 

266  plot,newc(0:2500) ,yrange*[-10. 0,10.0] ,y8tyle»l 

267  oplot,5.*longclean( 0:2500) ,color>120 

268  wait, 0.02 

269 

270  print, ' j«' ,strtrim(string( j) ,2) , '  k>' ,8trtrim( string (k) ,2) 
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271  andfor 

272 

273 

274  ;iirlijng*r«v«ra«(ahlft(20*alogl0(abs(ffsurf  )/magpinc) , 

275  ;ahifty,0));  tha  axptl  data 

276  nrlijfig«20*alogl0(aba(ffaurf )/iiagpinc)  ;  the  axptl  data 

277 

278  ;  tha  ahift  and  ravaraa  are  needed  to  make  my  NRL'a 

279  ;anglea  the  aame  aa  for  Benchmark  problema 

280 

281  bytx^bytacl (nrllmg,mln«-30 ,max»30 , top« I d. n_colora-l ) 

282 

283  waet,5 

284  tv,rebln(bytx,2*121,331),300,35 

285  xyouta, 0.4, 0.9 5, 'Experimental  data' , /normal 

286 

287  bytw-bytacl(20*alogl0(abn(modapec)/magpinc) ,minB~30,$ 

288  max*30,top«ld.n_colora-l) 

289 

290  wnet,5 

291  loadct,5 

292  tv,rebin(bytw,2*121,331),20,35 

293  xyouta, 0.08, 0.95, 'Windowed  model ', /normal 

294 

295  phaaec«atan( imaginary ( var ) ,float(var) ) 

296  phaaex«atan(imaglnary(ffaur£) ,float(ffaurf ) ) 

297  phaaew«atan(imaglnary(modapee) ,float(modapec) ) 

298  window, /free 

299  tvael,phaaew,50,50 

300  tvacl,phaaex,200,50 

301  tvael,phaaec,350,50 

302 

303 

304  goto,  jump 

305  hardcopy, f llename*atrtrim(  jobname, 2 )-•■$ 

306  atrtrlm(atrlng(  fix  (line) )  ,2)-l-'  .epa'  ,xaize>1.21,y8ize-3.31 

307  tv, byte 

308  cloaait 

309  hardcopy ,  filename" 'win '  -t-atnaidC  jobname ,1,2)+$ 

310  8trtrlm(8tring( fix (line) ) ,2)+' .epa' ,xaize"l .21,y8ize"3.3 

311  tv,rebln(bytw,121,330) 

312  cloaeit 

313  hardcopy , filename* ' ex ' tatrsiid ( jobnama ,1,2)+$ 

314  8trtrim(atrlng( fix (line) ) ,2)+' .epa' ,X8ize"1.21,y8ize"3.3 

315  tv,rebln(bytex,121,330) 

316  cloaeit 

317  jumpt 

318 

319  end 
1 
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